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ABSTRACT 

Supersonic  flow  past  oscillating  cylindrical  shells  and 
oscillating  flat-plate  cascades  is  analyzed  using  linearized 
characteristics  methods.   Pressure  distributions  and  gener- 
alized aerodynamic  forces  are  computed  for  arbitrary  radius 
to  length  ratios,  axial  and  circumferential  wave  numbers, 
Mach  numbers,  and  reduced  frequencies.   Excellent  agreement 
is  obtained  with  a  previous  solution  of  the  complete  un- 
steady linearized  potential  equation  by  Dowell  and  Widnall 
using  Laplace  transform  techniques.   Also,  pressure  dis- 
tributions are  computed  for  flat-plate  cascades  with  sub- 
sonic leading  edge  locus  for  arbitrary  solidity,  stagger 
angle,  frequency,  and  interblade  phase  angle.   For  compar- 
ison, a  two-blade  solution  is  developed  using  the  method 
of  singularities.   Good  agreement  is  obtained  for  this 
special  case. 
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I.    INTRODUCTION 

Supersonic  flutter  of  skin  panels  first  occurred  on  the 
German  V-2  missile  and  has  become  of  increasing  engineering 
importance  in  recent  years  with  the  development  of  super- 
sonic flight  vehicles.   Panel  flutter,  for  example,  was 
found  to  occur  on  the  X-15  Rocket  Aircraft  and  Saturn  Space 
Launch  Vehicle.   An  integral  part  of  the  flutter  analysis  of 
such  panels  is  the  determination  of  the  oscillatory  pressure 
distribution.   A  number  of  approaches  have  been  formulated 
which  resulted  mostly  in  approximation  theories.   A  detailed 
review  of  the  state  of  the  art  can  be  found  in  Reference  6. 
Only  one  "exact"  approach  was  worked  out  by  Doweil  and  Widnall 
using  the  complete  linearized  unsteady  potential  equation 
(Ref.  2).   However,  this  approach  requires  rather  complicated 
and  time-consuming  complex  integration  techniques. 

In  the  present  work  a  method  of  characteristics  procedure 
is  developed  to  determine  supersonic  flow  past  vibrating 
panels  and  cylindrical  shells  resulting  in  a  rapid  computation 
of  the  pressure  distribution  and  generalized  aerodynamic 
forces. 

A  second  flutter  problem  has  also  become  of  considerable 
interest  very  recently.   Due  to  the  continuing  trend  toward 
higher  tip  speeds  in  modern  jet  engine  fans  and  compressors, 
the  possibility  of  supersonic  blade  flutter  has  become  of 
considerable  concern.   A  recent  review  of  this  problem  area 
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has  been  given  in  "A  Workshop  on  Aeroelasticity  in  Turbo- 
machines"  (Ref.  13). 

In  the  present  study  the  method  of  characteristics  is 
applied  to  compute  supersonic  flow  past  a  finite  flat  plate 
cascade  with  subsonic  leading  edge  locus.   In  the  past,  cas- 
cade studies  have  proved  to  be  a  useful  tool  for  the  under- 
standing of  the  complicated  three-dimensional,  incompressible 
or  subsonic  flow  through  an  actual  compressor  stage.   Al- 
though it  is  doubtful  that  the  much  more  complex  mixed 
transonic-supersonic  flow  can  be  modeled  by  a  two-dimensional 
cascade  idealization,  the  analysis  of  an  infinitely  thin  flat 
plate  cascade  may  point  the  way  toward  a  fuller  understanding 
The  method  of  characteristics  again  proved  to  be  a  convenient 
tool  for  the  analysis  of  supersonic  flow  past  oscillating 
cascades . 
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II.    PROBLEM  FORMULATION 

Inviscid  adiabatic  flow  of  a  gas  along  the  positive  x- 
axis  is  considered  which  is  perturbed  by  either   A)  an  air- 
foil cascade  executing  small  pitch  oscillations  (Fig.  1)  or 
B)  small  amplitude  vibrations  of  the  surface  of  a  circular 
cylinder  whose  axis  is  aligned  with  the  free-stream  direction 
(Fig.  2). 

The  following  basic  equations  hold: 
The  equation  of  state  for  a  perfect  gas, 

p  =  PRT  (1) 

the  conservation  equation  for  mass, 

5i+  P(!i+  11+  *JL)  =  o  (2) 

Dt      3x    3y    3z 
the  conservation  of  momentum 

Du  +  1  8p_  =  Q 

Dt    p  3x 

Dv  +  1  3p  m    Q 
Dt    p  3y 

Dw  +  1  3_p  =  Q 
Dt    p  3z 

and  the  conservation  of  energy 

51=  o 

Dt 
D_ 

where  Dt  is  the  substantial  derivative  and  s  is  the  entropy 
of  the  fluid  particle. 

The  assumption  of  small  amplitude  vibrations  makes  it 
possible  to  introduce  the  small  perturbation  concept. 
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Figure  2.   Cylinder  Geometry 
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Therefore,  all  flow  quantities  are  assumed  as  small  pertur- 
bations linearly  superimposed  on  free-stream  quantities, 
i.e.  , 

u  =  u,,,  +  u  (5) 

V  =  v  (6) 

W  =  w  (7) 

5  =  0^  +  0  (8) 

and  the  pressure  and  density  perturbations  are: 

AP  =  P-Pco  (9) 

Ap  =  P-Pro  (10) 

The  pressure  perturbation  can  be  expressed  by  the  velocity 

of  sound  perturbation,  i.e., 

2 
P-Poo  =  rZi  PcOc0O(c-c0O)  (11) 

Specializing  to  planar  flow,  one  has  for  the  continuity 

equation: 

2   3c  .   2      9c       3u  ,     3v    n  ,,ox 

—  +  um  —  +  c^  —  +  cro  —  =0         (12) 


Y-l  at    y-1     3x      3x      3y 
and  the  Euler  equations: 

3u       3u     2      3c    A  ,,-,  v 

—  +  uro  —  +  Coo  —  =0  (13a) 

3t       3x    Y-l     3x 

3v  .     3v     2      3c    n  ,,  «,.» 

—  +  uro  —  +  cro  —  =0  (13b) 

3t       3x    y-1     3y 
Equations  (12)  and  (13)  are  the  governing  equations  for  the 
case  of  supersonic  flow  past  oscillating  cascades. 

For  the  case  of  supersonic  flow  past  vibrating  cylindrical 
shells,  the  equivalent  single  equation  for  the  perturbation 
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potential,  $(x,r,6,t)  in  cylindrical  coordinates  was  taken 

as  the  more  convenient  basic  governing  equation,  i.e., 

(1-M2)$   +$   +i$  +1  $  -2-  $  .--2*4.+.  =  0       (14) 
v      xx   rr  r  r  r*    69  c   xt  cz  tt 

GO  00 

For  a  derivation  of  this  equation  reference  is  made  to  one  of 
the  standard  texts,  e.g.,  Bisplinghof,  Ashley,  and  Half man, 
"Aeroelasticity" ,  (Ref.  1). 

For  a  complete  formulation,  proper  initial  and  boundary 
conditions  must  be  prescribed.   The  most  important  ones  are 
A)  the  flow  on  the  body  surface  must  be  tangent  at  any  instant 
in  time  and  B)  Sommerf eld's  radiation  condition  must  be  satis- 
fied, i.e.,  waves  must  propagate  away  from  sources  of 
disturbance. 

After  properly  linearizing  the  general  flow  tangency 
condition,  (see  Ref.  1)  one  obtains  for  the  flov;  tangency 
condition  on  a  thin  airfoil,  performing  a  time  dependent  small 
amplitude  motion  about  the  y=0  plane, 

dh        ah 

■m*t*"»:  at*=°  <15' 

where  h(x,t)  is  the  equation  of  the  airfoil  chord  line. 

For  the  vibrating  cylindrical  shell  the  flow  tangency 
condition  over  the  flexible  portion  becomes 

•J     -*+*!* 

r=R  3t    8x       O^!1 

*  I    =  0  x<0  (16) 

r=R 


15 


where  h(x,6,t)  in  this  case  describes  the  radial  shell 
deflection. 

With  the  assumption  of  a  perturbation  potential,  the 
momentum  equations,  equations  (3),  may  be  integrated  to 


+_2°_+u  +  — 5=c(t) 

3t   2     3X  ip 


(17) 


assuming  uniform  parallel  flow  far  upstream, 

c(t)=^» 
2 

Also,  the  density-pressure  relationship  can  be  expanded  as 


(18) 


p(p)=sp„+^)(p-p»)4^>= 


dp' 


(19) 


therefore, 


.P-Po 


•p  dp   =1 

"D   p  (p)   p   p   |_   p  Y_ 


OO      00 


dp 


~P-Pc 


(20) 


Hence,  the  perturbation  pressure  can  be  expressed  in  terms 
of  the  perturbation  potential,  i.e., 


P-P  =  ~P, 


3$  ,   3$ 
— +u  — 


(21) 


_9t    3x 
In  pressure  coefficient  form  this  becomes 

C  =P"Pqo  =  -2  $  -2$ 
P  I  "2    ^"2  t    U  X 

— p   U^         00        00 

^   00   00 


(22) 


For  further  analysis  it  is  convenient  to  use  non-dimen- 
sional quantities.   Introducing  for  the  oscillating  cascade 
case  the  blade  chord  as  the  unit  of  length  one  has  the 
following  non-dimensional  coordinates: 
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2  -  £,  y.-  £,  t  -  LJLi.  (23) 

c      c       c 

Furthermore,  introducing  the  assumption  of  harmonic  time- 
dependence,  one  can  define  a  reduced  frequency 

k  =  ^  (24) 

u 

where  w  is  the  circular  frequency  of  blade  oscillation. 

Following  previous  work  by  Teipel  (Ref.  3)  for  single 

blade  oscillations  the  non-dimensional  complex  amplitude 

functions 

.,,-  _,  iwt   u  /or-\ 

U(x,y)e     =  -  (25) 

TT,_  -v  iwt     1   v  n,, 

V(X,y)e     =  —  -  (26) 

_  ,-   -v  iwt     2   1c  /o-7\ 

C(x,y)e    -  —jjp  _  (27) 

are  introduced.   Substituting  these  equations  into  the  continu- 
ity  and  Euler  equations  (12)  and  (13)  gives 

|S  +  /m^T  |X  +  M2  |£  +  ikM2C  =  0     (28) 
3x  9y       3x 

|S  +  |£  +  ikU  =  0  (29) 

9x    3x 

~  +   -^-  |£  +  ikV  =  0  (30) 

where  the  overbar  on  x  and  y  is  being  dropped  for  convenience 
in  the  subsequent  analysis. 

Introducing  similar  non-dimensional  variables  with  dis- 
tances referred  to  flexible  cylinder  length,  velocities  to 
free-stream  speed  and  time  to  cylinder  length  divided  by 
free-stream  speed,  one  obtains  the  equation  for  the  non-dimen- 
sional perturbation  potential,  i.e., 
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(1-M2)$    +  $    +  -o    *nn     -  2M2$    -  M2*,.   =  0    (31) 
v      xx    rr    r*  69        xt      tt 

subject  to  the  boundary  condition 

$      ~   rr  +  t—  0<x<l 

r  r=R   3t    9X  (32) 

=  0  x<0 

If  the  radial  shell  deflection  is  assumed  as 

i  Tcf- 
h(x,6,t)  -  Z (x)  cos  ne-e  (33) 

where  k  is  now  referred  to  cylinder  length,  i.e., 

k  =  i^  (34) 

00 

Then  the  velocity  potential  must  also  have  the  form, 

ikt 
$(x,r,e,t)  =  <j>(x,r)  cos  n9»e         (35) 

Substitution  of  equation  (35)  into  equation  (31)  then  gives 
(1-M2U    +  f-  +  1*   -  Sy  +  -  2ikM2*   +  k2M2$  =  0   (36) 

AA      XX     xxx  X 

subject  to 

♦rl   p  =  M+ikZ         O^1  (37) 

r-R         =0  x<0 

The  pressure  coefficient  in  terms  of  the  non-dimensional 

velocity  potential  is 

C  =  -2[ik*  +  $  ]  (38) 

p  x 

and  the  pressure  coefficient  amplitude  is 

Cp  =  -2[ik<t>  +  <})  ]  (38a) 

These  are  the  basic  equations  used  in  the  subsequent  cylind- 
rical shell  vibration  analysis. 
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III.    THE  METHOD  OF  CHARACTERISTICS 

A.   CYLINDRICAL  SHELL 

In  Section  II,  equation  (36)  was  developed  as  the  partial 
differential  equation  for  <f>  (x,r)  that  requires  solution 
subject  to  the  boundary  condition,  equation  (37).   The  follow- 
ing developments  are  based  on  a  previous  paper  by  Platzer, 
Brix,  and  Webster,  (Ref .  7) . 

If  M>1,  equation  (30)  is  hyperbolic  and  has  real  character- 
istics, which  satisfy  the  ordinary  differential  equation 

(M2-l)dr2-dx2  =  0  (39) 

An  arbitrary  function  F(x,r)  has  the  following  derivatives 
along  the  two  characteristics: 


*L=F.H+F      dr=^EIF   +ip  (40) 

ds    x  ds     r  ds     M    x    M   r 

where  ds  is  the  differential  arc-length  along  the  character- 
istic ds2  =  M2dr2  and 


£l-  /S3  (4d 

ds     H  K      ' 

and 

—  =  +  -  (42) 

ds    "  M  *   ; 

Letting  Sj  be  the  arc- length  of  the  left-running  Mach  line 

and  s2  the  arc  -length  of  the  right-running  Mach  line,  one 

obtains  for  the  first  derivatives, 

Fi  =   ^r^  F  +  n  F  (43) 

1     M    x   M  r 
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and 

M   '  x  "  M  "  r 

and  for  the  second  derivatives 

-F 


Mz-1  1 

F,  =  "4t-=-  P._  -  sF. 


dF_L  =  I, 


12   '  ds2  '"  M2  [(l-M2)Fxx  +  Frr]  =  -F21 
Therefore,  equation  (30)  can  be  written  as 

d^T  "   2rM   /jpzj  (*i  +  *2)  +  (k   ^WU 


(44) 


(45) 


(46) 


Figure  3.   Cylindrical  Shell  Characteristics 

Network 


The  finite  difference  form  of  this  equation  for  the  points 
P  and  S  (see  Figure  3)  on  the  right-running  Mach  line  is 

4>1   (S)      -     jh    (P)     _     <j>i-<t>2     ~         .  t.      ,,      x       .       (V2        n2     .  M_, 

Similarly,  for  the  points  Q  and  S  on  the  left-running  Mach 
line  one  obtains 

r  Ui+*9)  +  (k2  — 

r2M2 


♦  ?(S)  -  »9(Q)  _  tl-^  _  JM      ■       2   n2 


(48) 
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Now  solving  for  the  values  of  <j>  \    and  ^  at  S  from  the  known 
values  at  P  and  Q  (Figure  3)  the  entire  flow  field  may  be 
computed  point  by  point.   For  this  purpose  the  values  for  <j>  j 
and  <J>2  in  equation  (47)  are  averaged  between  P  and  Sf  i.e., 

*x    =    *i(P)2+    MS)  (49) 

*2    ~        2~~^ (50) 


and 


<j)(S)    =    <j>(P)    +    1[<J)2(P)    +    <J>2(S)]AS  (51) 


is  obtained  by  integration  along  the  right-running  Mach  line. 

Similarly,  $  i    and  <j)2  in  equation  (48)  are  averaged  be- 
tween Q  and  S 

♦  ,  =  ^l(Q)2+  »1(S)  (52) 

♦  2  =  J2M    +   *?M         '  (53) 


and 


*(S)  =  <J>(Q)  +  I  t*i(Q)  +  *i(S)]  As        (54) 


2 

This  results  in  a  system  of  two  equations  for  ^(S)  and  <J>2(S), 

A*]  (S)  +  B(J>2  (S)  =  C  (55) 

D<J»!  (S)  +  E(|)2(S)  =  F  (56) 


21 


where 


n      As  ,         ikMAs  (57) 


4Mr(S)  2/M^l  4  M2r2(S) 

C   -    +l(P)[l   +    45^7^)         ^J=] 

+   #2(P)[.      **  _  _i^  +   Asi(k2    _  _n!_>] 

4Mr(P)  2/M^l  4  M2r2(S) 

(59) 

As  ikMAs         As2  ,,  o  n2         .  ,        ,,-.*.* 

D "  wnsr- "  j^5  *  "t- (k2-^ji^7)]    (60) 

E    =    -    1    -    j^y-    'JO**'  (61) 

4Mr(S)  2/M^I 

o   -    a     m.\  r    i     •       As                   ikMAs., 
F   =    (J>2  (Q)  [-1   -r   -        (G\    +   -— J 

4MrCQ;  2/M^T 

,_.  r  As  ikMAs         As2  n  o  n2         .  , 

4r(Q)M        2v^TT  4  M2r2(S) 

+  ♦(q)c-  ^(k2 si_)  -  ^.(k2  -  -a! — )i 

2  r2(S)M2  ^  M2r2(Q) 

(62) 
The  above  procedure  will  suffice  for  the  computation  of 
<j>  and  its  derivatives  at  a  general  flow  field  point,  but  needs 
to  be  modified  to  obtain  the  values  on  the  initial  Mach  line 
and  along  the  cylinder  surface. 

On  the  initial  Mach  line  (Figure  3)  it  is  known  that 
4>=<f>  i=0 .   Hence,  equation  (56)  becomes 

E<fr2(S)  =  F  (63) 
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where  E  is  known  from  equation  (61)  and  F  simplifies  to 

F=  *a(Q)  t-1  +  35EW  +  p=5=l  (64) 

On  the  cylinder  surface  the  normal  velocity  is  prescribed, 
i.e.  , 

4>r  -  |  Ui-M  =  ff  +  ikz  <65> 

whereas  equation  (55)  must  hold  along  the  right-running  Mach 
line.   The  values  for  <$>  i  and  ^  on  the  surface  are  then  found 
from 

c  + 1%  (M +  ikz)B 

HIS)    =  L^ (66) 

and 

C  -  1  (§£  +  ikZ)B 

*2(s)  =  Vt-b (67) 

The  pressure  distribution  over  the  vibrating  cylindrical 
shell  is  obtained  from  the  linearized  pressure  amplitude, 
equation  (3  8a) ,  and  the  generalized  aerodynamic  force  Q   is 


found  from 


where 


0   =  1   f   C   (x)  V  (x)dx  (68) 

ror   2  o  Pm     r 


Z(x)  =  Z  Am  Y  (x)  (69) 

m  m 
m 

and  C   (x)  is  the  aerodynamic  pressure  amplitude  due  to  the 

deflection  ¥  (x) . 
m 

These  equations  were  programmed  in  Fortran  IV.   The 
program  description  appears  in  Section  IV-A. 
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B.   CASCADE 

In  Section  n  equations  (28)  -  (30)  were  developed  as  the 
governing  system  of  partial  differential  equations  for  the 
analysis  of  the  cascade  problem.   A  solution  of  these  equations 
can  be  obtained  using  the  method  of  characteristics.   The 
following  development  of  this  method  is  largely  based  on  the 
previous  work  by  Teipel  (Ref.  3)  and  Chalkley  (Ref.  9),  where 
further  details  can  be  found. 

Transforming  the  system  into  a  new  coordinate  system, 
£  (x,y)  and  n(x,y),  and  considering  now  a  system  of  three  simul- 
taneous equations  in  unknowns  U  ,  V  ,  and  C  ,  the  equations 


may  be  written 


2-  n  tt  „..  "  cot  a  V_n. 


%    U   +  cot  a  5  V,  +  MZC  Cr  =  -  U  ny    ~~   ~   n  Y 

-  M2C  n   -  ikM2C  (70) 

n  x 

CUr  +  £C=-Un   -  C  n   -  ikU         (71) 
x  K  x  £      nx    nx 

K   V  +  tan  a  £  C,  =  -  V  n   -  tan  a  C  n   -  ikV   (72) 
x  £         y  £      nx         n  y 

If  U  ,  V  ,  and  C   are  made  indeterminate  across  lines  of 

r\   =  constant,  the  determinant  of  the  matrix  of  coefficients 

of  the  above  equations  must  equal  zero,  or 

-  C3  -  I     U  2  -  M2  C2]  =  0  (73) 

x    x   y        x 

The  solutions  to  this  equation  give  the  three  character- 
istic directions  in  the  physical  (x,y)  plane: 

P-  =    0  (74) 

dx 


^  =  ±  tan  a  (75) 
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In  order  to  find  the  compatibility  relations  of  the  first 
derivatives  of  U,  V,  and  C  along  the  characteristics,  equations 
(70)  -  (72)  are  written  along  lines  of  5  -  constant/   n  = 
constant,  and  along  streamlines,  that  is,  for 

str(x,y)  =  constant,  -r%-   =  0        (76) 

3  X 

£(x,y)  =  constant,  -£  =  (77) 

9x  /m^T 

8  v    -  1 
ri(x,y)  =  constant,  ~-  =  (78) 

If  an  arbitrary  function  of  x  and  y  in  the  equations  is  de- 

3  f 
noted  by  f  then,  - —  holding  str  constant  becomes 

o  X 

'II' str  =  i  (79) 

3f 

Likewise,  - —  holding  £  constant  becomes, 

3  X 

(ii)  =  i£  +  -i —  ii  (80) 

3f 
and,  - —  holding  n  constant  becomes 

a  X 

3x  n    3x   /M2_i  3Y 
Combining  these  equations  and  rewriting, 

£  =    <M>str  .  (82) 

H  =  i[(i>5+  <!fV  (83' 

||  =  i/iP=l    [<H),    -    (|f)  J  (84) 

dy        2  °x  s  3x   n 

Substituting   these   relations   into   the   equations   of  motion 
gives 
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1      r,3U.  ,3U,     ,  1      ,„2     i  >      ,,3V,  ,DV,     , 

+  iM2[<ll's  +    'f'n1    +  ik>l2c  =  °  (85) 

I   I(f'?+    <fV    +f   [(I§>C    +    'I'n1    +iku=°         <86> 

[<i'5+    <fV    +T    "I'e    "    <IV    +   ^=0  (87) 

Subtracting  equation  (86)  from  equation  (85)  gives  an  equation 
that  is  first  added  to  equation  (87)  and  then  subtracted  from 
equation  (87) .   The  result  is  two  compatibility  relations  that 
must  be  satisfied  along  lines  of  constant  £  and  lines  of 
constant  ri,  i.e., 

^K    +    <1§>*  +  ik[V  +  ■J—  (M2C-U)]  =  0         (88) 
3x  5     3x  *  M2-l 

(^n  "  <!§>„  +  ik[V  "  ^—  (M2C-U)]  =  0         (89) 
dX  n    c3x  n         M2-l 

The  third  compatibility  relation  is  obtained  by  considering 
equation  (86)  along  a  streamline 

(|~)  .   +  (~)  .   +  ikU  =  0  (90) 

3x  str     3x  str 

If  irrotational  flow  is  assumed,  the  equation  of  irrota- 

tionality 

2H  -  H.  =  0  (91) 

3y    3x  v   ; 

may  be  substituted  for  the  second  equation  of  motion  in  the 
system.  '  With  the  assumption  of  simple  harmonic  motion,  it 
may  be  written  in  terms  of  Teipel's  amplitude  functions, 

|£  -v^^I^  =  0  (92) 

3y        3x 

The  non-dimensional  system  of  equations  can  then  be  written 

as, 
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|£  +/  Mz_x  3V    M2  3£    ikM2c  =  (93) 

ax        ay      ax 

|£  +  |£  +  ikU  =  0  (94) 

ax   ax 

9£  -/  m^-1  |^  =  0  (95) 

ay       ax 
By  transforming  equations  (93)  -  (95)  to  the  arbitrary  co- 
ordinate system,  %    and  r\ ,    the  system  can  be  considered  one  of 
simultaneous  equations  in  U  ,  V  ,  and  C  .   By  making  any  of 
these  derivatives  indeterminate,  the  results  of  setting  the 
determinant  of  the  matrix  of  their  coefficients  equal  to  zero 
are  equations  (74)  and  (75),  the  characteristic  directions. 

By  writing  the  equations  of  motion,  as  before,  along 
lines  of  %    -   constant  and  r\   =  constant  and  performing  the 
same  algebraic  manipulations,  the  two  compatibility  relations 
can  be  determined: 

(!£),  ~  (|2)   +  ik  -Hi-  (U-C)  =  0         (96) 

dX  t,       3X  t,         M2-l 

(|£)   +  (|X)   +  ik  -Hi-  (U-C)  =  0         (97) 

ax'n   vax  n     m2_1 
The  third  compatibility  relation  is,  as  before, 

(~)  4.   +  (~)  ,_   +  ikU  =  0  (90) 

3x  str    ax  str 

The  equations  of  motion  are  now  reduced  to  a  system 
containing  only  derivatives  with  respect  to  x.   It  is,  thus, 
possible  to  solve  these  equations  using  finite  differences  in 
the  following  manner:   separating  the  real  from  the  imaginary 
parts,  the  equations  can  be  written, 
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9UR         dCR 
( — -)     +  (t—^)     -  kUT  =  0 
9x   str     9x  'str     I 


(98) 


9U  9  C 

^str  +  <35T>str  +  kUR=  ° 


(99) 


3U, 


3  V, 


(jrr)r         i^r),    -   k  JSi-  (U   -  C  )  =  0 

9x.  £     9x   £      M2-l 


(100) 


3U        3V  2 

^5  -  (^5  +  kfr  (ur-cr>  =  °     (101) 

M.  —  1 

3UR        8VR  M2 

(~)   +  (r=S)n  -  k  -^—  (U_  -  C  )  =  0    (102) 

dx   n    9x   n     m2-1 

1 

3U        3V  2 

(^)n  +  (r-i)   +  k  -2—  <Up  -  C.)  -  0    (103) 

dX   f|      dX   r)       M2_i     K      K 

To  write  these  equations  in  finite  difference  form  the 
computational  molecule  shown  in  Figure  4  will  be  used.   With 
the  velocity  values  all  known  at  P,,/  P]2/  an<3  p21'  a^  t^ie 
values  at  P~?  can  be  calculated.   The  distance  Ax  is  arbitrary. 


P22 


-%* 


Figure  4.   General  Flow  Field  Molecule 
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If  F  is  an  arbitrary  flow  quantity,  its  partial  deriv- 
atives can  be  written  as 

F   -  F 
(2L)  22    xll 

l3x  str     Ax  (104) 

F    -  F 

(|I)   =  21 li  (105) 

8X  5       J-  AX 

C||)n  -  ^  "  ?21  (106) 

9X  n       [AX 

The  values  of  the  flow  quantities  in  the  equations  are 
averages ,  such  that 

(F)str  "  \    <F11  +  F22>  <107> 

(F)5  =7  (F12  +  F22'  (108) 

(F)„  "  \    (F21  +  F22'  (109' 
With  these  relations  substituted  into  equations  (98)  -  (103) , 
the  system  becomes  after  some  algebraic  manipulation, 

U22R  +  C22R  "  AIU22I  =  K12R  -  (110) 

U22I  +  C22I  +  AIU22R  =  K12l  (111)  - 

U22R  "  V22R  +  BI(C22I  "  U22l}  =  K34R  (112) 

U22I  "  V22I  +  BI(U22R  "  C22R)  =  K34I  (113) 

U22R  +  V22R  +  BI(C22I  "  U22I}  =  K56R  (114) 

U22I  +  V22I  +  BI(U22R  "  C22R)  =  K56I  (115) 


where 


K12R   U11R  +  CHR  +  AiUlll  (116) 

K12I  =  U11I  +  C1H  "  AiUllR  (117) 
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K34R=  U12R~  V12R+  BI (U12I  "  C12I>  (118) 

K34I    U12I  "  V12I  *"  BI(U12R  "  C12R)  (119) 

K56R  :  U21R  +  V21R  +  BI (U21I  "  C21I)  (120) 

K56I    U21I  +  V21I  "  BI(U21R  "  C21R)  (121) 


and 


A_  =  j   kAx  (122) 

1      M2 
B=[k  — Ax  (123) 

1    k         M2-l 
As  shown  in  Teipel  (Ref.  3),  solving  these  equations 
gives, 

U22R  -  (l-AIBI)[^K34R+K56R)-DIK12I]+2BllI(K34I+K56I)+BIK12R] 

(1-Aj.B^.)2  +  UBj)2  (124) 

U22I   =    (1-Vl)tr'K34I+K56I)+BIK12R1-2BI[T'K34R+K56R)-BIK12I1 

(1-AjBj)2    +     (2BJ.)2  (125) 

V22R  -  7  (K56R  "  K34R>  <126> 

V22I  "  7  <K56I  "  K34l'  <127) 

C22R  "  K12R  "  U22R  +  AIU22I  (128) 

C22I  =  K12I  "  U22I  "  flIU22R  (129) 

These  are  the  finite  difference  equations  for  a  point  in  the 
general  flow  field. 

At  a  point  on  the  upper  surface  of  a  cascade  blade,  the 
normal  velocity  V"22  is  prescribed  by  the  movement  of  the  blade 
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This  velocity  was  given  in  the  problem  formulation.   Here  the 
computational  molecule  is  shown  in  Figure  5. 


y  A 


Figure  5.   Upper  Surface  Molecule 
The  applicable  equations  are, 

U22R  +  C22R  ~  AIU22I    K12R 


U22I  +  C22I  +  AIU22R  "  K12I 


U22R  "  BIU22I  +  BIC22I  "  K56R  "  K34R 

U22I  +  BIU22R  "  BIC22R  =  K56I  "  K34l 

where  K.  ~  and  Kcr    are  as  before,  but, 
1Z  bo 


K34R  :  V22R 


K34I    V22I 


(130) 
(131) 
(132) 
(133) 


(134) 
(135) 


as  given  by  the  boundary  conditions.   Solving  for  the  flow 
quantities  gives 


U22R  =  (1-AIBI)[K56R-K34R-BIK12I]+2BI[K56I-K34I+BIK12R] 


(1-AjBj)2  +  (2BI)2 


(136) 


U22I    *    (1-AIBI)[K56rK34I+BIK12R^2Bi^56R-K34R-BIK12I] 


(1-AjBj)2    +     (2BX)2 


(137) 


C22R    K12R  "  U22R  +  AIU22I 


(138) 
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and  V22R  and  V22J 


C22I  =  K12I  "  U22I  "  AIU22R         (139) 
are  known  from  the  boundary  conditions. 
These  are  the  finite  difference  equations  for  the  flow 
velocities  on  the  upper  surface  of  a  cascade  blade. 

Likewise,  at  the  lower  surface  of  a  cascade  blade,  the 
normal  velocity,  V^  is  prescribed  by  the  movement  of  the 
blade,  and  was  given  in  the  problem  formulation.   The  compu- 
tational molecule,  different  from  that  at  the  upper  surface 
of  the  blade,  is  shown  in  Figure  6. 

y  11 


-**— 


Figure  6.   Lower  Surface  Molecule 
Here,  the  applicable  equations  are, 


U22R  +  C22R  '  AIU22I    K12R 


U22I  +  C22I  +  AIU22R  "  K12I 


U22R  "  BIU22I  +  BIC22I  ~  K56R  +  K34R 


U22I  +  BIU22R  "  BIC22R  "  K56l  +  K34I 


where  K12  and  K  .  are  as  originally  stated,  but, 


(140) 
(141) 
(142) 
(143) 


K56R  :  V22R 


K56I    V22I 


(144) 
(145) 
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as  given  by  the  boundary  conditions.   Thus,  the  finite 
difference  equations  for  the  flow  velocities  at  the  lower 
surface  of  a  cascade  blade  are, 

0     _  d-A^)  IIWK34R-BIK12I]+2BI[K56I+K34I+BIK12R1 
22R  (1-AjBj.)2  +  (2BX)2  (146) 

U    =  (1-AIBI)[K56I+K34I+BIK12R]-2BI[K56R+K34R-BIK12I1 
221  (1-Aj.Bj)2  +  (2BJ)2  (147) 

C22R  =  K12R  "  U22R  +  AIU22I  (148) 

C22I  "  K12I  "  U22I  '  AIU22R  (149) 

where  V?2   and  V-_   are  given  from  the  boundary  condition. 

To  obtain  the  conditions  on  the  initial  left-running 
Mach  line,  assume  that  P.. ,  and  P.--,  as  shown  in  Figure  4  are 
just  in  the  free-stream  and  then  let  Ax  shrink  to  zero. 
Since  K..  2  =  K  fi  =  0  and  A   and  B_  go  to  zero  as  Ax  goes  to 
zero,  the  initial  finite  difference  equations  for  £=constant 
are , 

U22  +  C22  =  0  (150) 

U22  "  V22  +  V12  "  °  -.  (151) 

U22  +  V22  =  0  (152) 


thus, 


U22  -  "  V22  "  "  C22  (153) 


Equation  (96)  can  be  written, 

9xJ  % 


(|^)r  =  -  ikU  (154) 
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Integrating  gives, 


M2 

-ikT  ~x 


U  =  Ul   ne  X^M2-1  (155) 

1  x=0 

thus , 

M2  M2 

U„D  =  -Voop(0)cos  (k  ~ —  x)-Vo9T(0)sin  (k  -ii —  x)  (156) 
22R      22R  M2-l  M2-l 

U22I  =  ""V22I(0)COS  (k    x)+V22R(0)sin  (k    x)  (157) 

M2-l  M2-l 

V22R  "  -  U22R  (158) 

(159) 

(160) 

(161) 

For  the  initial  right-running  Mach  line,  P. ,  and  P..  ~  of 
Figure  4  are  assumed  to  be  just  in  the  free-stream  and  Ax 
is  brought  to  zero.   Since  K,2  =  K-..  =  0,  and  A  =  E  =  0, 

U22  +  C22  =  0  (162) 

U22  -  V22  =  0  (163) 

U22  "  U21  +  V22  "  °  (164) 


V22I 

=  "  U22I 

C22R 

=  ~  U22R 

c 

^221 

=  "  U22I 

thus, 


U22  "  V22  "  "  C22  (165) 


Equation  (97)  can  then  be  written, 


8x  n 
Integrating  gives, 


(|^)   =  -  ikU  (166) 


■  v   M2 
U  =  V(0)e~1K-^-[x  (167) 
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Thus, 

M2  M2 

U„D  =  V_o_(0)cos  (k  — —  x)+V99T(0)sin  (k  — —  x)   (168) 
22R     22R  M2-l  M2-l 

0  0 

U0,T  =  V,,T(0)cos  (k  — — x)-V.,„(0)sin  (k  — — x)   (169) 
221    221  M2-l  M2-l 

V22R  "  U22R  (170) 

C22R  "  "U22R  <171> 

V22I  "  U22I  <172> 

C22I  =  "U22I  <173> 

At  the  first  position  encountered  on  a  blade  other  than 
the  first  blade  there  are  two  grid  points  for  which  the  flow 
field  quantities  must  be  computed;  one  on  the  upper  surface 
and  one  on  the  lower  surface  of  the  blade . 

As  in  Figure  7  the  blade  is  assumed  to  protrude  a  very- 
small  distance  past  P^/  thus  necessitating  the  computation 
of  the  upper  and  lower  flow  velocities  at  this  point.   In  this 
case,  the  equations  developed  earlier  for  an  upper  surface 
grid  point  are  used  for  P^pn  an<^   t^lose  developed  for  a  lower 


surface  grid  point  are  used  at  POOT 


Blade 
Figure  7.   New  Blade  Molecule 


35 


At  a  position  just  aft  of  a  blade  there  are  also  two 
grid  points  to  be  accounted  for;  one  above  and  one  below  the 
slip  plane  in  the  wake.   As  in  Figure  8  the  point  P??  is 
assumed  to  be  a  very  small  distance  aft  of  the  blade's  trail- 
ing edge. 


Wake 


Figure  8.   Wake  Molecule 


At  this  position  the  equations  for  an  upper  blade  grid 
point  apply  to  PppTI  and  those  for  a  lower  grid  point  apply  to 


221/ 


i.e.  , 


TT  +    C  -     A    IT 

22RU  22RU  I    22IU 


K 


12RU 


U22IU  +  C22IU  +  AIU22RU    K12IU 


U22RU  +  BI(C22IU  "  U22IU}    K56R  "  V22RU 


U22IU  +  BI(U22RU  ~  C22P 


^   —  T<f      —  V 

561     22IU 


(174) 
(175) 
(176) 
(177) 


and 


U22RL  +  C22RL  "  AI°22IL  =  K12RL 


U22IL  +  C22IL  +  AIU22RL  "  K12IL 


U22RL  +  BI(C22IL  "  U22IL)  =  V22RL 


(178) 

(179) 

+  K34R  (180) 


U22IL  +  BI(U22RL  "  C22RL)  =  V22IL  +  K34I  (181) 
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where 

K12RU  -  U11RU  +  C11RU  +  VllIU  (182> 

K12IU  =  U11IU  +  C11IU  "  AIU11RU  (183) 

K12RL  =  U11RL  +  C11RL  +  VllIL  U84) 

K12IL  =  "llIL  +  C11IL  "  AIU11RL  <185> 

K34R=  U12R  "  V12R+  BI(U12I  "  C12I>  <186> 

K34I  "  U12I  "  V12I  '  VU12R  "  C12r'  (187' 

K56R  =  U21R  +  V2iR  +  VU21I  "  C21I>  (188) 

K56I  "  U21I  +  V21I  "  BI(U21R  -  C21R>  (189) 

Equations  (174)  -  (181)  form  a  system  of  eight  equations  with 
twelve  unknowns.   However,  since  it  is  known  that  the  wake 
cannot  support  a  pressure  difference  and  that  the  normal 
velocities  on  either  side  of  the  wake  must  be  equivalent,  i.e., 

V22RU  "  V22RL  ^  V22R  (190) 

V      =  V      =  V  (191) 

22IU     22IL     221  K         } 

C  =  C      =  C  (192} 

^22RU    ^22RL  "  ^22R  K^    } 

C      =  C      =  C  (193) 

^22IU    ^22IL  "  ^221  v    } 

we  can  eliminate  four  of  the  unknowns  and  achieve  a  determinate 

system  of  equations.   Equations  (174)  -  (181)  become 

U22RU  +  C22R  "  AIU22IU  =  K12RU  (194) 

U22IU  +  C22I  +  AIU22RU  =  K12IU  ("5) 

U22RU  +  V22R  +  BI(C22I  "  ^HJ*  "  K56R       (196) 

U22IU  +  V22I  +  BI(U22RU  "  C22RJ  =  K56I       (197) 
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U22RL  +  C22R  "  AIU22IL  =  K12RL  (198) 

U22IL  +  C22I  +  V22RL  "  K12IL  (199) 

U22RL  "  V22R  +  VC22I  "  U22IL»  =  K34R       (200) 

U22IL  "  V22I  +  VU22RL  "  C22R>  =  K34I       (201) 

The  pressure  distributions  along  the  surfaces  of  the  blades 
can  be  determined  from  equation  (11) .   Dividing  both  sides  of 
this  equation  by  p   and  using  Teipel's  amplitude  function  for 
c,  the  equation  can  be  written, 

P(x,y)  =  YM2C(x,y)|y=0-^  yssQ-  (202) 

where 

P(x,y)eiwt  =  E^~  (203) 

Thus  the  non-dimensional  pressure  distributions  along  the 
upper  and  lower  surfaces  of  the  blades  in  a  cascade  can  be 
computed  in  a  point  by  point  calculation  of  U,  V,  and  C  along 
lines  of  £=constant  and  n-constant. 

These  equations  were  programmed  in  Fortran  IV.   The  pro- 
gram description  appears  in  Section  IV-B. 

C.   THE  LINEARIZED  METHOD  OF  CHARACTERISTICS  OF  SAUER  AND  HEINZ 

Sauer  and  Heinz  (Ref.  11)  developed  a  linearized  method 
of  characteristics  for  supersonic  flow  past  bodies  of  revo- 
lution at  zero  angle  of  attack.   This  method  also  is  derived 
from  the  linearized  potential  equation  for  steady  supersonic 
flow,  i.e.  , 

(M2-lHxx  -  4>rr  -  I  *r  =  0  (204) 
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Sauer  and  Heinz  showed  that  the  compatibility  relations  can 
be  written  as 

2£tf-'-/H*=rr  |£  (205) 

d  t,  d  t, 


and 


where 


AlEXi   =   /m^T  r  |H  (206) 

3  x]  3  n 


u  -  1^  (207) 

d  X 


v  =  |±  (208) 

d  r 


and 


5  s  x  -  /M2-l  r  (209) 


n  =  x  +  /M2-l  r  (210) 

are  the  characteristic  coordinates  as  explained  previously 

Writing  this  equation  in  finite  difference  form  one 
obtains  for  u  and  v  at  point  4  (Figure  9)  using  the  known 
values  in  points  2  and  3 ,  - 


(rv)  ->-  (rv)  0+/M2-l  (r  u0+rru0 
U4  =  3     2- H-J L_£  (2ii) 

/M2-l  (rr+r  ) 


rr(rv)-+r  (rv) 0-/M2-l  rrr  (u.-u.) 
(vr)   =  -i 2 2 £ L-2 1 — L.         (212) 

r   +  r 
5    n 

where  r   and  r   are  the  average  r  coordinates  between  3  and 
4  and  between  2  and  4  respectively.   Computing  the  flow  pro- 
perties at  the  body  using  the  exact  flow  tangency  condition 
leads  to 
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/M2-l  r  u.+  (rv).-r,  (^-)  _ 

n  4      4   5   dx  5  ,on  -jv 

u   =   (213) 

D 

dh. 


/m2-1  rn  +  r5(57>5 


v5  "  (§>5(1+u5>  (214> 

where  r  =  h(x)  describes  the  body  contour. 

However,  in  order  to  apply  the  Sauer-Heinz  procedure  to 
the  case  of  an  infinitely  long  cylinder  with  a  wavy  surface 
between  0<x<l  a  different  starting  procedure  is  needed.   It 

is  easily  seen  from  Figure  9  that 

h        1 
r.?    =  ^0  ( +  1)  (215) 

2    /Mz-1  tan  6 


realizing  that 


where 


and 


R   +  R 

r   =  -i~ 216) 

n      2 


R1  =  RQ  +  r±  (217) 


R2  =  RQ  +  hQ  (218) 


one  obtains 

ho        1 

rn  =  R  +  TT      (3  + )  <219) 

n        4        v/MT^T  tan  6Q 

Since  there  is  no  disturbance  upstream  "of  the  leading  edge 
Mach  line,  one  obtains  from  equation  (213) 

-  (R+hn)  tan  6n 
U2  =   — — (220) 

/M2-l  r  +(R+h.)  tan  6A 

n     o      o 

where  the  point  2  is  the  first  body  point. 
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Finally,  the  pressure  coefficient  is  obtained  from 

C   -  2u  -  v2  +  u2  (M2-l)  (221) 

P 

The  program  listing  for  these  equations  is  given  in 
Appendix  I. 


Figure  9.   Sauer-Heinz  Characteristic  Network 
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IV.    COMPUTATIONAL  PROCEDURE 

A.   CYLINDRICAL  SHELL  COMPUTER  PROGRAM 

This  program  calculates  the  pressure  distribution  over  a 
harmonically  vibrating  cylindrical  shell.   It  is  valid  only 
for  a  cylinder  in  a  supersonic  flow  field  vibrating  with 
small  amplitude  vibrations. 

The  program  calculates  the  flow  field  properties  adjacent 
to  the  cylinder  using  the  method  of  characteristics  finite 
difference  equations  developed  previously  in  Section  III. 
The  pressure  distribution  is  then  integrated  using  the  trape- 
zoidal rule  according  to  equation  (68)  to  obtain  the  general- 
ized aerodynamic  force. 

The  computational  diagram  shown  in  Figure  3  is  used  to 

calculate  the  flow  field  quantities  at  S  in  the  general  flow 

field  and  on  the  cylinder  surface.   The  distance  As  shown  in 

the  figure  is  determined  by 

L  cos  a  /ooo\ 

As  =  (222) 

v 

where  v  is  the  grid  fineness  ratio.   The  grid  fineness  ratio 
is  an  arbitrary  input  variable  into  the  program  and  is  equal 
to  one  less  than  the  total  number  of  grid  points  on  the 
cylinder  surface. 

Input  parameters  are  entered  into  the  program  on  two  data 
cards.   The  first  card  contains  the  date,  e.g./  11  APR  1973, 
in  the  first  twelve  columns.   The  second  card  contains  the 
following  input  data: 
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SPACE  DATA 

1-20  Free-stream  Mach  number,  M 

21-30  Radius  to  length  ratio,  R 

31-40  Reduced  frequency,  k 

41-45  Axial  mode  number,  m 

46-50  r,  as  defined  in  equation  (68) 

for  the  generalized  aero- 
dynamic force. 

51-55  Circumferential  mode  number,  n 

56-60  Grid  fineness  ratio,  v 

where  the  first  three  quantities  are  floating  point  numbers 
and  the  rest  are  integers.   The  definitions  of  these  para- 
meters are  given  again  in  the  alphabetical  listing  of  all 
program  variables  in  Appendix  A. 

The  program  has  five  subroutines.   They  are  as  follows: 

COMPXR  Computes  the  x  and  r  position 

of  the  given  grid  point. 

MACHLN  Computes  the  flow  quantities 

along  the  initial  left-running 
Mach  line  at  the  given  grid 
point. 

GENFPT  Computes  the  flow  quantities  at 

a  general  flow  field  grid  point. 

PANPNT  Computes  the  flow  quantities  on 

the  cylinder  surface  at  the  given 
grid  point. 

INTEG  Computes  the  generalized  aero- 

dynamic force. 

Flow  quantities  along  the  initial  left-running  Mach  line 

are  computed  using  the  equations  derived  in  the  method  of 

characteristics  for  initial  conditions.   Those  quantities  on 

the  cylinder  surface  are  computed  using  the  equations  derived 
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for  the  boundary  condition.   The  values  of  the  velocity 
potential  and  its  derivatives  are  derived  as  a  result  of  the 
analysis  of  the  flow  tangency  condition  and  have  been  given 
previously.   The  procedure  starts  at  the  initial  Mach  line 
and  continues  down  successive  right-running  Mach  lines.   The 
pressure  coefficient  is  computed  when  the  cylinder  surface 
is  met  and  stored  in  an  array  for  later  use  in  the  integration 
of  the  generalized  aerodynamic  force. 

A  flow  diagram  of  this  program  is  given  in  Appendix  B. 

B.   OSCILLATING  CASCADE  PROGRAM 

This  program  calculates  the  pressure  distribution  over 
the  upper  and  lower  surfaces  of  the  blades  of  a  subsonic 
leading  edge  locus  cascade  which  is  oscillating  in  pitch  at 
arbitrary  frequency  and  interblade  phase  angle.   It  is  valid 
only  for  a  cascade  in  a  supersonic  flow  field  vibrating  with 
small  amplitude  vibrations.   The  program  calculates  the  flow 
field  quantities  within  the  area  governed  by  the  aft-running 
Mach  lines  emanating  from  the  leading  edge  of  the  first  blade 
and  the  forward-running  Mach  lines  emanating  from  the  trailing 
edge  of  the  last  blade  of  the  cascade  (see  Figure  10) .   These 
flow  field  properties  are  calculated  using  the  method  of 
characteristics  finite  difference  equations  developed  previously 
in  Section  III-B. 

The  computational  molecules  shown  in  Figures  4,  5,  6,  7 , 
and  8  are  used  to  calculate  the  flow  field  properties  in  the 
general  flow  field,  upper  surfaces,  lower  surfaces,  leading 
edge  of  a  staggered  blade,  and  the  wake  aft  of  a  blade, 
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respectively.   The  distance  Ax  shown  in  the  figures  is 
determined  by 


2d 

Ax  =  —  cot  a 


(223) 


where  v  is  the  grid  fineness  ratio.   This  parameter  is  an 
arbitrary  input  variable  into  the  program  and  is  equal  to  one 
less  than  the  total  number  of  grid  points  on  any  Mach  line 
that  runs  from  one  blade  to  the  next. 


Figure  10.   Cascade  Computational  Area 

Input  parameters  are  entered  into  the  program  by  means  of 
the  NAMELIST  option  of  Fortran.   In  using  this  option,  input 
parameters  are  combined  in  a  titled  list  and  then  referred 
to  in  the  program  by  this  title  (NAM1) .   The  real  advantage 
of  using  NAMELIST  is  that  each  data  card  need  only  contain 
the  value  of  one  input  variable  while  the  format  of  the  data 
card  is  such  that  the  variable's  name  is  punched  on  the  card 
as  well  as  its  value.   In  this  way,  when  a  different  value 
for  an  input  variable  is  desired,  rather  than  changing  the 
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entire  data  deck,  only  one  easily  identifiable  data  card  need 
be  changed. 

The  format  for  a  NAMEL.IST  data  deck  is  as  follows:  The 
first  column  in  each  card  of  the  data  deck  is  left  blank. 
The  first  card  starts  in  the  second  column  with  the  symbol  & 
(ampersand)  followed  immediately  by  the  title  of  the  NAMELIST 
(here,  NAM1) .   The  next  cards  list  the  variables  and  their 
input  values,  one  per  card,  in  any  order.   The  format  here, 
starting  in  the  second  column  and  followed  by  a  comma,  is: 
variable  name  =  value.   The  last  card  in  the  data  deck  and  the 
signal  to  the  computer  that  the  NAMELIST  has  ended  again 
starts  in  the  second  column  with  the  symbol  &  immediately 
followed  by  END. 

The  value  of  an  input  variable  may  be  written  in  any 
format  as  long  as  integer  variables  are  written  without  deci- 
mal points  and  real  variables  are  written  with  a  decimal  point. 
Blanks  are  taken  as  zeroes,  but  a  comma  must  appear  somewhere 
between  a  desired  input  value  and  the  next  data  card.   The 
date  of  the  run  is  entered  on  a  data  card  in  the  first  twelve 
spaces.   This  data  card  is  inserted  prior  to  the  NAMELIST 
dataset. 

To  illustrate  the  input  procedure,  the  following  example 

is  proposed: 

v  =  100        d  =  0.8 

M  =  /2  3  =  55° 

k  =  0.75        6  =  45° 

xQ  =  0.5 

5  blades 
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The  data  cards  needed  are  as  follows: 


CARD 

FORMAT 

1 

22  JUNE  1973 

2 

&NAM1 

3 

NGRDFN=100, 

4 
5 
6 
7 


FSTRMN=1. 414214, 
REDFRQ=0.80, 
XSUBO=0.50, 
TNWDST=0.80, 

STGANG=55.0, 


DESCRIPTION 

Date 

NAMELIST  title 

Grid  fineness  ratio, 
must  be  an  even  integer, 
subject  to  the  constraint 
of  equation  (227)  . 

Free-stream  Mach  number 

Reduced  frequency 

Elastic  axis  position 

Vertical  distance  between 
blades . 

Complimentary  stagger  angle, 
tan  3  >  cot  a . 

Interblade  phase  angle 

Number  of  blades 

End  of  NAMELIST 

The  definitions  of  these  parameters  are  given  again  in  the 

alphabetical  listing  of  all  program  variables  in  Appendix  C. 

The  program  has  nine  subroutines.   They  are  as  follows: 

INPUT      Reads  in  all  input  data  and  sets  up  the  finite 
difference  grid. 

INTIAL     Initializes  all  the  flow  quantities  at  the  leading 
edge  of  the  first  blade  and  initializes  most  of 
the  logic  variables. 

MACHLN     Computes  the  flow  quantities  along  the  initial  left 
and  right-running  Mach  lines  at  the  given  grid 
point. 

TOP        Computes  the  flow  quantities  on  the  upper  surface 
of  a  blade  at  the  given  grid  point. 


9 

FAZE=45.0, 

10 

NUMBLD=5 , 

11 

SEND 
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GENFPT     Computes  the  flow  quantities  at  the  given  general 
flow  field  point. 

BOTTOM     Computes  the  flow  quantities  on  the  lower  surface 
of  a  blade  at  the  given  grid  point. 

NEWBLD     Computes  the  flow  quantities  at  the  leading  edge 
of  all  staggered  blades. 

WAKE       Computes  the  flow  quantities  in  the  wake  aft  of 
each  blade  at  the  given  grid  point.   This  sub- 
routine uses  IBM  subroutine  SIMQ  for  the  solution 
of  its  simultaneous  equations. 

COMPXY     Computes  the  x  and  y  position  of  the  given  grid 
point. 

In  INPUT  the  grid  is  determined  using  the  grid  fineness 
ratio,  blade  spacing,  number  of  blades,  and  the  free-stream 
Mach  number.   The  value  of  Ax  as  shown  in  equation  (223)  is 
obtained.   The  stagger  angle  is  then  adjusted  such  that  the 
distance  the  upper  blades  are  staggered  back  is  the  nearest 
multiple  number  of  Ax  increments.   The  "compatible  stagger 
angle"  is  then  printed  out.   The  blades  may  be  staggered  such 
that  the  first  staggered  blade  is  at  least  one  Ax  increment 
aft  of  the  initial  left-running  Mach  line,  but  no  further 
forward ,  or, 

d  •  (tan  6  -  cot  a)>Ax  (224) 

Flow  quantities  along  the  initial  left  and  right-running 
Mach  lines  are  computed  using  the  equations  derived  in  the 
method  of  characteristics,  Section  B,  for  initial  conditions 
Those  flow  quantities  along  the  blades'  surfaces  are  computed 
using  the  equations  derived  for  the  appropriate  boundary 
condition.   The  equations  for  the  normal  velocities  are  de- 
rived as  a  result  of  the  analysis  of  the  flow  tangency  con- 
dition and  have  been  given  previously.   The  flow  quantities 
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in  the  wake  aft  of  each  blade  are  computed  using  the  set  of 
simultaneous  equations  derived  as  a  result  of  the  assumptions 
made  about  the  flow  quantities  on  either  side  of  the  wake  and 
are  given  in  the  Method  of  C haracteristics ,  Section  III-B. 

The  procedure  starts  at  the  leading  edge  of  the  first 
blade  and  computes  the  flow  quantities  at  grid  points  along 
the  initial  right-running  Mach  line.   It  then  jumps  back  to 
the  initial  left-running  Mach  line  and  computes  the  flow 
quantities  along  the  next  right-running  Mach  line,  using  the 
applicable  subroutines  at  the  appropriate  grid  points.   Only 
flow  quantities  on  the  Mach  line  presently  being  computed 
and  the  Mach  line  which  preceeded  it  are  kept  in  storage  at 
any  one  time.   Flow  quantities  at  successive  grid  points  are 
computed  in  like  manner  until  the  trailing  edge  of  the  last 
blade  of  the  cascade  is  met. 

A  flow  diagram  of  this  program  is  given  in  Appendix  D. 

C.   SAUER-HEINZ  COMPUTER  PROGRAM 

This  program  computes  the  pressure  distribution  over  a 
cylindrical  shell  or  body  of  revolution.   It  is  written  so 
that  the  shape  of  the  body  or  cylindrical  surface  may  be 
easily  changed.   It  is  valid  only  for  non-oscillating  bodies 
or  shells  in  a  supersonic  flow  field  at  zero  angle  of  attack. 
The  following  discussion  applies  to  the  cylindrical  shell. 

The  program  calculates  the  flow  field  adjacent  the 
cylinder  using  the  method  of  characteristics  finite  difference 
equations  developed  previously  in  Section  III-C. 
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The  computational  diagram  shown  in  Figure  9  is  used  to 

compute  the  flow  field  properties  at  4  in  the  general  flow 

field  and  at  5  on  the  cylinder  surface.   The  distance  As 

shown  in  the  figure  is  determined  by, 

,     cos  a  /ooc\ 

As  =  (225) 

where  v  is  the  grid  fineness  ratio.   The  grid  fineness  ratio 
is  an  arbitrary  input  variable  into  the  program  and  is  equal 
to  one  less  than  the  number  of  grid  points  on  the  initial 
Mach  line.   The  length  As  ,  n=l,2,3...  is  determined  by  the 
intersection  of  the  right-running  Mach  lines  with  the  surface 
Input  parameters  are  entered  into  the  program  on  one  data 
card.   This  card  contains  the  following  input  data: 
SPACES  DATA 

1-10  Free-stream  Mach  number,  M 

11-20  Maximum  amplitude  of  the  axial  mode,  hn 

21-30  Cylinder  radius  to  length  ratio,  R. 

31-40  Axial  mode  number,  m 

41-50  Grid  fineness  ratio,  v 

where  the  first  three  quantities  are  floating  point  numbers 

and  the  last  two  are  integers.   The  definitions  of  these 

quantities  are  given  again  in  the  alphabetical  listing  of 

all  program  variables  in  Appendix  E. 

The  program  has  three  subroutines.   They  are  as  follows: 

COMPXR  Computes  the  x  and  r  position  of  the 

grid  point  and  determines  the  position 
of  the  intersection  of  the  current 
right-running  Mach  line  with  the  body. 

GENFPT  Computes  the  flow  quantities  at  a  general 

flow  field  grid  point. 
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SURFAC  Computes  the  flow  quantities  at  the 

given  surface  point. 

Flow  quantities  on  the  initial  left-running  Mach  line 
are  assumed  to  be  zero.   Those  quantities  on  the  surface  are 
computed  using  the  appropriate  equations  developed  in 
Section  III-C,  for  surface  grid  points. 

The  procedure  starts  at  the  initial  Mach  line  and  con- 
tinues down  successive  right-running  Mach  lines.   The  pressure 
coefficient  is  computed  when  the  surface  is  met. 

A  flow  diagram  of  this  program  is  given  in  Appendix  F. 
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V.    RESULTS  AND  COMPARISONS 

A.   CYLINDRICAL  SHELL 

Pressure  distributions  and  generalized  aerodynamic  forces 
were  computed  for  sinusoidal  deflection  modes,  i.e., 

*  (x)  =  sin(mTTx)  m=l,2,...       (226) 

Figures  11  through  16  show  results  for  the  generalized  aero- 
dynamic forces  0, ,  and  Q, ~  as  a  function  of  circumferential 
mode  number  for  two  radius  to  length  ratios,  r=0.5  and  r=1.0, 
at  a  Mach  number  of  /2~  and  reduced  frequencies  of  zero  and 
one.   Also  shown- is  the  comparison  with  two  previous  solutions 
by  Dowell  and  Widnall  (Ref .  2)  and  by  Matsuzaki  and  Kobayashi 
(Ref.  5).   The  Dowell  and  Widnall  solution  is  also  based  on 
the  complete  linearized  unsteady  potential  equation,  but  uses 
Laplace  transform  techniques  for  its  solution.   As  can  be  seen, 
agreement  between  this  solution  and  the  method  of  character- 
istics solution  is  excellent.   For  the  sake  of  interest,  a 
recent  approximate  solution  by  Matsuzaki  and  Kobayashi,  valid 
only  for  large  values  of  the  circumferential  mode  number  is 
also  shown. 

Figures  17  through  24  show  the  generalized  aerodynamic 
forces  Q-,-|R  and  Q,  ~      as  a  function  of  radius  to  length  ratio 
again  for  a  Mach  number  of  /2,  zero  reduced  frequency,  and 
circumferential  mode  numbers  n=0,l,2,  and  3.   The  main  purpose 
of  these  figures  is  to  show  that  the  infinite  cylinder  theory 
is  a  poor  approximation  for  small  radius  to  length  ratios,  but 

approaches  the  exact  theory  for  large  values  of  R . 
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Figure  18.   Generalized  Aerodynamic  Force  vs 

Radius  to  Length  Ratio 
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Figure  22.   Generalized  Aerodynamic  Force  vs. 

Radius  to  Length  Ratio 


3.0 


64 


♦ 


1.6 


1.4 


1.2 


Q 


1.0 
0.8 


Infinite  Cylinder  (Ref .  4) 

Q  Dowel 1  and  Widnall  (Ref.  2) 
Present  Study 


M=/2 

n=2 

k=0 


0.6 


0.4 


0.2 


0.0 


/ 


-»-*- 


0.0 


Figure  23 


1.0 


R 


2.0 


Generalized  Aerodynamic  Force  vs 
Radius  to  Length  Ratio 


3.0 


65 


4 


1.8 


1.6 


1.4 


1.2 


1.0 


Q 


12 


R 
0.8 


0.6 


0.4 


0.2 


0.0 


• Infinite  Cylinder  (Ref.  4) 

□  Dowell  and  Widnall  (Ref.  2) 
Present  Study 


M=/2 

n=3 

k=0 


0.0 


1.0 


2.0 


3.0 


Figure  24.   Generalized  Aerodynamic  Force  vs. 

Radius  to  Length  Ratio 


66 


Finally,  Figure  25  shows  the  change  in  pressure  distri- 
bution as  a  function  of  radius  to  length  ratio  for  the  case 
M=/2~,  m=l,  n=0,  k=0.   These  distributions  were  computed  using 
the  characteristics  method  outlined  in  Section  III-A  and  the 
Sauer-Heinz  theory  as  described  in  Section  III-C.   Complete 
agreement  was  found  between  the  two  methods. 

B .   CASCADE 

Figures  26  through  29  show  the  instantaneous  pressure 
distributions  on  the  first  and  second  blades  of  a  finite 
cascade  for  four  different  values  of  reduced  frequencies. 
The  pressure  difference  between  lower  and  upper  blade  surface 
is  plotted  for  the  instant  when  the  airfoil  is  at  its  mean 
horizontal  position  and  pitching  nose  up.   Results  for  these 
cases  were  given  very  recently  by  L.  E.  Snyder,  of  Pratt  and 
Whitney  Aircraft  (Ref.  8),  for  the  following  parameter 
combinations : 

M  =  1.5 

Solidity  =  1.33  (d  =  0.317) 

Stagger  angle  =  25°  (B  =  65°) 

xQ  =  0.5 

S  =  40°   (Blade  2  leads  blade  1) 

k  =  0.04,0.332,0.832,1.25 
It  can  be  seen  that  appreciable  differences  exist  in  a 
number  of  cases,  in  particular  in  the  wave  reflection  regions 
These  differences  are  not  fully  understood  at  this  time. 
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Figure  26.   Cascade  Pressure  Coefficient  Difference 
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Snyder  uses  a  finite  difference  procedure  and  further  work  is 
required  to  evaluate  the  validity  and  accuracy  of  this  pro- 
cedure is  compared  to  the  method  of  characteristics  procedure 
developed  in  the  present  work.   However,  a  check  could  be  made 
for  the  first  blade  solutions  for  which  the  method  of  singular- 
ities of  Garrick  and  Rubinow  (Ref.  12)  is  available.   This 
solution  was  programmed  according  to  the  equations  given  in 
Appendix  H  and  the  results  are  shown  in  Figures  26  through 
29.   The  method  of  characteristics  solution  is  in  complete 
agreement  with  the  Garrick-Rubinow  prediction  whereas  consider- 
able deviations  are  seen  to  exist  for  the  finite  difference 
predictions . 

A  further  check  using  the  method  of  singularities  was 
developed  for  two  blades  oscillating  with  zero  phase  difference 
This  approach  is  also  described  in  Appendix  H  and  a  comparison 
with  the  method  of  characteristics  is  given  in  Table  1  below: 

TABLE  1 


k 

0 

.1 

0 

.5 

M 

/2 

1 

.8 

/2 

1 

.8 

X 

M.O.C. 

M.O.S. 

M.O.C. 

M.O.S. 

M.O.C.    M 

.O.S. 

M.O.C. 

M.O.S. 

0.1 

.979 

.979 

.667 

.660 

.693 

.694 

.640 

.486 

0.2 

.977 

.977 

.666 

.661 

.686 

.687 

.638 

.511 

0.3 

.974 

.974 

.666 

.662 

.680 

.681 

.635 

.536 

0.4 

.972 

.972 

.666 

.663 

.676 

.677 

.631 

.562 

0.5 

.969 

.969 

.666 

.664 

.674 

.675 

.627 

.588 

0.6 

.966 

.966 

.665 

.666 

.675 

.676 

.623 

.614 

0.7 

.963 

.963 

.665 

.667 

.678 

.679 

.618 

.641 

0.8 

.959 

.959 

.665 

.668 

.682 

.683 

.614 

.668 

0.9 

.956 

.956 

.664 

.669 

.689 

.690 

.607 

.694 
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Further  method  of  characteristics  results  for  a  finite 
subsonic  leading  edge  cascade  are  shown  in  Figures  30  through 
32  for  the  following  parameter  combinations: 

M=/2         xQ=0.5 

d=1.0         k=0.5 

0=59.8°        6=120° 
As  mentioned  previously,  the  computer  program  is  written 
sufficiently  general  so  that  the  pressure  distribution  can 
be  computed  and  printed  out  for  arbitrary  specified  Mach 
number,  reduced  frequency,  stagger  angle,  blade  spacing,  and 
number  of  blades  subject  to  the  following  constraint: 

v  tan  a  [ (N-l)  (tan  6-cot  a)  +  i]<800     (227) 

Equation  (227)  is  the  expression  which  governs  the  maximum 
grid  fineness  that  may  be  used  with  a  given  set  of  input 
parameters.   This  constraint  is  necessary  because  the  program 
is  dimensioned  such  that  the  maximum  number  of  grid  points 
allowed  on  any  one  right-running  Mach  line  is  400. 
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VI.    CONCLUSIONS  AND  RECOMMENDATIONS 

A  solution  of  the  complete  linearized  unsteady  potential 
equation  has  been  developed  for  supersonic  flov;  past  vibrating 
cylindrical  shells  using  the  method  of  characteristics. 
Pressure  distributions  and  generalized  aerodynamic  forces 
were  computed  for  shells  of  arbitrary  radius-to-length  ratios; 
axial  and  circumferential  wave  numbers,  Mach  number,  and  re- 
duced frequencies.   A  solution  for  this  same  problem  was  pre- 
viously obtained  by  Dowell  and  Widnall  using  Laplace  transform 
techniques.   Excellent  agreement  was  obtained  in  all  cases. 
It  is  believed  that  the  present  method  can  be  applied  more 
easily  than  methods  requiring  complex  variable  techniques. 
It  should  be  noted  that  the  procedure  can  easily  be  modified 
to  include  arbitrary  axial  deflection  modes. 

A  further  solution  of  the  complete  lineariz.ed  unsteady 
potential  equation  using  the  method  of  characteristics  was 
developed  for  supersonic  flow  past  vibrating  flat  plate 
cascades  with  subsonic  leading  edge  locus.   Pressure  distri- 
butions were  computed  for  finite  cascades  of  arbitrary 
stagger  and  interblade  phase  angles  performing  small  ampli- 
tude pitch  oscillations.   Preliminary  comparisons  with  a 
recent  finite  difference  solution  by  Verdon  and  Snyder  showed 
reasonable  agreement.   However,  further  work  is  necessary  to 
assess  the  validity  and  accuracy  of  the  present  method.   This 
should  be  done  by  further  comparisons  with  the  Verdon  and 
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Snyder  procedure  and  Fourier  and  Laplace  transform  methods 
currently  under  development  at  Stevens  Institute  of  Technology, 
General  Motors  Corporation  (Detroit  Diesel,  Allison  Division) , 
and  General  Electric.   Also,  flutter  calculations  should  be 
carried  out  using  the  method  of  characteristics  predictions 
of  the  oscillatory  aerodynamic  forces  which  would  allow  com- 
parisons with  the  Pratt  and  Whitney  cascade  flutter  tests. 
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APPENDIX  A 
CYLINDRICAL  SHELL  COMPUTER  PROGRAM  VARIABLES 

A.  LOGIC  VARIABLES 

The  following  variables  appear  in  the  logic  statements 

of  the  program: 

IHAVEP         The  number  of  the  computation  along  any 

right-running  Mach  line.   IHAVEP  is  used  as 
the  first  subscript  for  a  grid  point. 

ISWTCH         Switch  variable  used  to  designate  the 

right-running  Mach  line  being  calculated, 
alternately  set  equal  to  "1"  or  "2"  and 
used  as  a  second  subscript. 

JSWTCH         Switch  variable  used  to  designate  the 

previously  computed  right-running  Mach 
line,  always  the  opposite  of  ISWTCH  and 
used  as  a  second  subscript. 

KOUNT  The  value  of  IHAVEP  at  the  cylinder  surface. 

ML  The  number  of  grid  points  along  the  initial 

Mach  line  and  along  the  cylinder  surface. 

B.  QUANTITY  VARIABLES 

The  following  variables  in  the  program  take  on  the  values 

defined  below.   All  are  floating  point  quantities  unless 

otherwise  stated.   The  dimensioned  variables  are  listed  first. 

AMAGCP(400)       The  imaginary  part  of  the  pressure  coefficient 

CP(400)  The  complex  pressure  coefficient. 

PHI (400,2)        The  derivative  of  the  velocity  potential 

along  the  right-running  Mach  line. 

PHK(400,2)       The  derivative  of  the  velocity  potential 

along  the  left-running  Mach  line. 

QPHI  (400,2)       The  velocity  potential. 

80 


REALCP (400,2) 
X,R(400,2) 
A,B,C,D,E,F 
BETA 

DELTAS 
DSTSTR 
FINGRD 

HDSTRL 

HEIGHT 

II 

IKM 

IKW 

IPIII1,    IPHI2 

K 

MACH 

MANGL 

MFREQ 

NFREQ 

NR 

PARW 

PI 
PX 

RADIUS 
W 


The  real  part  of  the  pressure  coefficient. 

The  x  and  r  coordinates  of  a  grid  point. 

See  equations  (57)  -  (62) 

The  square  root  of  the  quantity  Mach  number 
squared  minus  one. 

The  step  size  along  the  characteristics. 

The  horizontal  distance  between  grid  points. 

The  grid  fineness  ratio.   An  even  integer 
less  than  400. 


One-half  DSTSTR. 


The  increment  HDSTRL  ?   tan  a . 

The  square  root  of  minus  one. 

M 
The  factor,  lk 

M2-l 

The  quantity,  ik  sin  kittx 

The  velocity  potential  as  described  in 
equations  (51)  and  (54) . 

The  reduced  frequency. 

The  free-stream  Mach  number 

The  Mach  angle,  a. 

The  axial  mode  number.   An  integer. 

The  circumferential  mode  number.   An  integer 

r,  as  defined  in  equation  (68).   An  integer. 

The  partial  derivative  of  the  deflection 
taken  with  respect  to  x. 

3.141593 

The  partial  derivative  of  the  velocity 
potential  taken  with  respect  to  x. 

The  radius  to  length  ratio. 

The  cylinder  deflection. 
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APPENDIX  B 


CYLINDRICAL  SHELL  FLOW  DIAGRAM 


Read/Write  DATE 

Read  Input  Data 

Define  XLENGT,  HEIGHT,  DELTAS, 
DSTSTR,  HDSTRL,  ML, 
T  RNGLH,  MANG L,  II,  BETA 


ISWTCH=  2 
JSWTCH=  1 
IHAVEP=  1 
KOUNT=  1 


Compute  QPHI(1,1),  PHK(1,1),  X(l,l), 
R(l,l),  PI,  PARW,  PHI (1,1), 
CP ( 1 , 1 ) ,  REALCP ( 1 ) ,  AMAGCP ( 1 ) 

Write  Input  Data 

KOUNT=  2 


104 

CALL  COMPXR 

Compute  PARW,  W,  XKM,  I KM,  XKW,  IKW 


(iHAVEP.EQ.l)- 


ICALL  MACHLN 


(lHAVEP  .  EQ  .  KOUNt}- 


1 


CALL    PANPNT 


(CALL    GENFPT| 


50 

Compute  PX,  PR 


Q 


1000 

IHAVEP.EQ.KOUNT 


> 


25 


IHAVEP=  IHAVEP  +  l] 


8.2 


f25  \ 


ISWTCH=1 
JSWTCH-2 
IHAVEP=1 


K0UNT==K0UNT+1| 


((KOUNT-1)  .EQ.MLJp- 


ISWTCII=2 
JSWTCH=1 
IHAVEP=1 

_J 
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APPENDIX  C 


OSCILLATING  CASCADE  COMPUTER  PROGRAM  VARIABLES 


A.   LOGICAL  VARIABLES 

The  following  variables  appear  in  the  logic  statements  of 
the  program: 


ICO 


I HAVE P 


ISWTCH 


JLINE 


JSWTCH 


LCOUNT 


LIMIT 


LIMW 


MAXI 


NEWDST 


Set  equal  to  "1"  when  a  blade  or  wake  com- 
putation is  made.   Set  equal  to  "0"  after 
the  next  general  flow  field  computation 
occurs . 

The  number  of  the  computation  along  any 
right-running  Mach  line  starting  with  zero 
at  the  initial  left-running  Mach  line.   The 
value  IHAVEP+1  is  the  first  subscript  for 
a  flow  quantity. 

Switch  variable  used  to  designate  the 
right-running  Mach  line  being  calculated, 
alternately  set  equal  to  "1"  or  "2"  and 
used  as  a  second  subscript. 

The  Mach  line  counter  incremented  each  time 
a  right-running  Mach  line  is  completed. 
Reset  equal  to  the  value  of  IHAVEP  at  the 
leading  edge  of  a  new  blade. 

Switch  variable  used  to  designate  the 
previously  computed  right-running  Mach 
line,  always  the  opposite  of  ISWTCH  and 
used  as  a  second  subscript. 

The  value  of  IHAVEP  at  a  surface  or  wake 
grid  point. 

The  maximum  value  of  IHAVEP. 

The  value  of  LCOUNT  at  the  trailing  edge 
of  a  blade. 

The  value  of  JLINE  on  a  right-running  Mach 
line  which  intersects  a  blade's  leading 
edge. 

The  number  of  Ax  increments  in  a  blade  chord 
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NSTPTS 

HUM 

NUMBLD 
NUMOLD 
NUMPI 
OLDL 


The  number  of  Ax  increments  between  the 
initial  Mach  line  and  the  leading  edge 
of  the  first  blade. 

Set  equal  to  one  less  than  the  number  of 
the  blade  encountered  on  a  Mach  line. 

The  number  of  blades  in  the  cascade. 

Preserves  the  maximum  value  of  NUM. 

NUMOLD  plus  one. 

Preserves  the  value  of  LCOUNT  at  the  highest 
blade  or  wake  yet  encountered. 


B.   QUANTITY  VARIABLES 

The  following  variables  in  the  program  take  on  the  values 

defined  below.   All  are  floating  point  quantities  unless 

otherwise  stated.   The  dimensioned  variables  are  listed  first. 

A  (8,  8)  Matrix  of  coefficients  in  wake  computation. 

B(8)  Right-hand  side  of  equations  (194)  -  (201) 

U22R,U22l,V22R,V22l,C22R,C22I (400,2) .   The  velocities,  real 

and  imaginary  of  general  flow  field,  upper 
blade,  and  upper  wake  points.   Equivalent 
to  the  respective  Teipel  amplitude  functions, 

U22R  throu9h  C22I* 

U33R,U33I,V33R,V33I,C33R,C33I (400, 2) .   Same  as  above,  but 

used  only  at  lower  blade  and  lower  wake 
points . 


X,Y(400,2) 

AI 

ARG 

BI 

DELTA 

DELTAS 

DSTSTR 


The  x  and  y  coordinates  based  on  chord  length 
Ay  as  defined  in  equation  (122) . 


The  argument  k 


M- 


M' 


B   as  defined  in  equation  (123). 

The  interblade  phase  angle  in  radians. 

Incremental  distance  between  grid  points 
along  Mach  line. 


The  increment  Ax. 
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FAZE 
FMANGL 
FNGDPT 
FSTRMN 
Gl ,  .  .  .  , G9 

HDSTRL 

K12R, . . . ,K56I 

KS 

NGRDFN 

REDFRQ 
STAG 

STGR 

STGANG 

S 

TNWDST 

TRNGLH 

U 


VIPANL 


VRPANL 


The  interblade  phase  angle  in  degrees. 

Mach  angle,  a. 

Floating  point  NGRDFN. 

The  free-stream  Mach  number. 

Intermediate  factors  used  in  velocity 
computations  in  GENFPT,  TOP,  BOTTOM, 
and  NEWBLD. 

The  increment  Ax/2. 


1?R,...,KrfiT  as  defined  in  equations  (116)  - 


K 

(1 


Dummy  variable  used  in  SIMQ.   An  integer. 

The  grid  fineness  ratio.   An  even  integer, 
less  than  400. 

The  reduced  frequency. 

The  distance  the  second  blade  is  staggered 
back  from  the  y-axis. 

The  distance  from  the  initial  Mach  line  to 
the  second  blade's  leading  edge. 

The  stagger  angle  in  degrees. 

The  factor  /M*-l. 

The  vertical  distance  between  adjacent  blades 

The  increment  yAx  cos  a . 


The  factor  cos (k 


The  factor  sin(k 


M- 


M2-l 

M2 
M2-l 


x)  . 


x). 


W 


The  imaginary  part  of  the  normal  flow 
velocity  at  the  first  blade's  leading  edge. 

The  real  part  of  the  normal  flow  velocity 
at  the  first  blade's  leading  edge. 

M2 


The  factor 


/M^T 
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XSUBO  The  elastic  axis  position  based  on  chord 

length. 

XLNGTH  The  length  of  the  initial  Mach  line  between 

the  planes  of  adjacent  blades. 

XN  Floating  point  NUM. 

XNEW  The  x-coordinate  on  a  blade  based  on  the 

blade's  leading  edge. 
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APPENDIX  D 


CASCADE  FLOW  DIAGRAM 


Read/Write  DATE 

CALL  INPUT 

OLDL=0 

ICO=0 

CALL  INTIAL 

Compute  NEWDST,  LIMIT, 

LIMW,  MAXI 
NUM=0 
NUMOLD=0 
lNUMPI=l 


1000 

CALL  COMPXY 


(lcount  .  eq  .  0} 


('iHAVEP  .  EQ  .  LCOUNT  > 


t  r 


(iHAVEP.EQ.o)-^- 


IHAVEP.EQ.NUMPI*NSTPTS 
.AND. 
JLINE.EQ.MAXI 


CALL  GENFPT 


(lHAVEP  .  EQ  .  LIMIT)- 


5 

t  > 


IHAVEP . EQ . IHAVEP+1 


T 


[CALL    MACHLNl 


(lHAVEP  .  EQ  .  LIMIt)^— [T 


|lHAVEP  =  IHAVEP  +  H 


CALL    NEWBLD 
ICO=l 

NUMOLD=NUM 

LC0UNT=IHAVEP4NGRDFN 

JLINE=IHAVEP 

IHAVEP=IHAVEP+1 

OLDL= IHAVEP 

NUMPI=NUM+1 

MAXI=NGRDFN+NUMPI *NSTPTS 

NUM=NUM-1 


2 
ICO=l 


(J 


IHAVEP . GE . LIMW 


yt. 


CALL  WAKE 


T 
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CALL    TOP 
CALL    BOTTOM 


(iHAVEP  .  EQ  .  LIMIT)-3- |ICO  =  0 


(iHAVEP  .  EQ  .  LIMIT*)— — 


NUM.GE.NUM3LD-1 
.AND. 
IIAVEP.EQ.  LIMIT 


IHAVEP . EQ . IHAVEP+H 


(iHAVEP  .  EQ  .  LIMIt)-^ 


(num.eq.o) 


F 


NUM=NUM-1 

LCOUNT=LCOUNT+NGRDFN 

LIMW=LIMW-NSTPTS 


|IHAVEP  =  fiIAVEP  +  l| 


(num.eq.o) 


T 


LCOUNT=LCOUNT+NGRDFN 
NUM=NUM - 1  


IHAVEP=0 

NUM=NUMOLD 
JLINE=JLINE+1 


(iSWTCH.EQ.l) 


T 


F 


ISWTCH=1 

JSWTCH=2 


ISWTCH=2 
JSWTCH=1 


■«<- 


(numold.gt.o) 


T 


LCOUNT^OLDL 

LIMW=NUMOLD*NSTPTS+NEWDST 

OLDL=OLDL+l 


JLCOUNT . EQ . LCOUNT+1 


STOP 
END 
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APPENDIX  E 

SAUER-HEINZ  COMPUTER  PROGRAM  VARIABLES 

A.  LOGIC  VARIABLES 

The  following  variables  appear  in  the  logic  statements 

of  the  program: 

IHAVEP  The  number  of  the  computation  along  any 

right-running  Mach  line  starting  with  zero 
at  the  initial  left-running  Mach  line.   The 
value  IHAVEP+1  is  the  first  subscript  for  a 
flow  quantity. 

ISWTCH  Switch  variable  used  to  designate  the 

right-running  Mach  line  being  calculated, 
alternately  set  equal  to  "1"  or  "2"  and 
used  as  a  second  subscript. 

JSWTCH  Switch  variable  used  to  designate  the 

previously  computed  right-running  Mach 
line,  always  the  opposite  of  ISWTCH  and 
used  as  a  second  subscript. 

KOUNT  The  value  of  IHAVEP  at  the  surface. 

B.  QUANTITY  VARIABLES 

The  following  variables  in  the  program  take  on  the  values 

defined  below.   All  are  floating  point  quantities  unless 

otherwise  stated.   The  dimensioned  variables  are  listed  first. 

CP(200)  The  pressure  coefficient. 

DELX,DELR(200)    The  incremental  distances  between  grid 

points  on  a  right-running  characteristic. 

U  (200,2)  The  streamwise  velocity  perturbation. 

V(200,2)  The  product  vr . 

X,R(200,2)  The  x  and  r  coordinates  of  the  grid  point. 

A,B  Used  as  factors  in  the  equations  for  u  and  v. 
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BETA 
DELTAS 

DELTAX 

DELTAY 

FINGRD 

HO 

MACH 

MANGL 

MFREQ 

NGFN 

NUMPTS 

PI 

RETA 

RO 

RZETA 

S 

SLP 

XM 

XMAX 

YO 


The  square  root  of  the  quantity  M2-l. 

The  incremental  distance  between  adjacent 
grid  points  along  the  initial  Mach  line. 

The  change  in  x  from  one  point  to  the 
next  on  the  initial  Mach  line. 

The  change  in  r  from  one  point  to  the  next 
on  the  initial  Mach  line. 

Floating  point  NGFN. 

The  maximum  amplitude  of  the  deflection. 

Free-stream  Mach  number. 

The  Mach  angle,  a. 

The  axial  mode  number,  an  integer. 

The  grid  fineness  ratio,  an  even  integer, 
less  than  200. 

The  number  of  grid  points  on  the  initial 
Mach  line. 

3.141593 


The  radius  to  length  ratio. 

The  length  of  the  initial  Mach  line. 

Slope  of  the  surface  at  a  surface  grid  point 

Floating  point  MFREQ. 

The  value  of  x  at  a  surface  grid  point. 

The  intercept  of  a  right-running  Mach  line 
with  the  r-axis. 
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APPENDIX  F 


SAUER-HEINZ  FLOW  DIAGRAM 


Read/Write  Input  Data 
Compute  BETA,  MANGL,  FINGRD, 
NUMPTS,  S,  DELTAS, 
DELTAX,  DELTAY,  XM 


ISWTCH=2 
JSWTCH=1 
IHAVEP=0 
I=IHAVEP+1 


Compute  X(l,l),  R(l,l),  U(l,l),  V(l,l) 
KOUNT=l    


CALL    COI4PXR 


(iHAVEP.NE.p}^- 


JT\T 


(lHAVEP  .  EQ  .  KOUNt) 


U(1,ISWTCH)=0 
V(1,ISWTCH)=0 
IHAVEP=IHAVEP+1 


CALL  SURFAC 


CALL  GENFPT 
IHAVEP=IHAVEI 


+  lj  (x.GE.l.o)^- 


IIHAVEP^SI 


STOP 
END 


JSWTCH  .  EQ  .  lJ-= 

F 

ISWTCH=11 
JSWTCH=2 

ISWTCH=2 
JSWTCH=1 

KOUNT=KOUNT+ll 


i-«- 


92 


APPENDIX  G 

SUPERSONIC  FLOW  PAST  AN  INFINITELY  LONG 
WAVY  WALL  CYLINDER 

Leonard  and  Hedgepeth  (Ref .  4)  gave  a  solution  for  super- 
sonic flow  past  an  infinitely  long  wavy  wall  cylinder.   De- 
tailed computations  using  this  solution  were  carried  out  by 
Anderson  (Ref.  14).   For  the  cylinder  wall  deflection 

h  =  hftsin  mux  cos  n  6  (G-l) 

where  m  is  the  number  of  half-waves  per  unit  axial  length 
and  n  indicates  the  number  of  circumferential  waves,  one 
obtains  for  the  pressure  distribution. 

miT 


Cp   =    2hn   ~  A,     sin(mTvx    +    y.  )  (G-2) 

r  0  M        1  1 


where 


] 


iy,               M         H       2     [mTiR(M2-l)  2] 
A,  e      1    =      — — —     n 

1  /Mz-1  i_ 

H      2'  [mTTR(M2-l)  2]  (G-3) 

n 

For    the    actual    computation   equation    (G-3)    was    rewritten   as 

A  ei,rl   =  Jn(Y)    -   i  Yn(Y) ^ 

1  'Jn-l(^-Jn+l^>]-i[Vl(^-Yn+l^)]       M^ 


where 


(G-4) 


Y    =   irnrR    /Mz-1  (G-5) 

For    the   case   n   =    0    the    equation   was    rewritten   as 


m  i»i  -  jq(y)  : L  yq(y)      M 

A   e      1        •      — (G-6) 

-Jx(y)    +   i  yjl(y)  Jm^i 


93 


and 


Al    =    lAlelVFll  (G_7) 


¥      =    tan"1    Im(AielH'1) 

Yl  an       [ — rj, ]  (G-8) 

Re(A1e1^l) 


The  program  listing  appears  in  Appendix  I 
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APPENDIX  H 


TWO  BLADE  SOLUTION  USING  THE 
METHOD  OF  SINGULARITIES 


As  shown  by  Garrick  and  Rubinow  (Ref.  12),  the  velocity 
potential  of  a  single  blade  oscillating  in  supersonic  flow  is 


M' 


/  \  ~1         ,X-y/Mz-l         ,,_»      -ikrrj — t-(x-C)t     r     ^M      ft 7T2 7772 — T~\ T-i   if 

d)  (x,y)  =       /  J  v(£)e    M^r     J«  I /(x-U  ^~  (M/~D  y  ]  d? 

^F^T  0  U  M2-l 


(H-l) 


-*-d/M2~l->^    / 


d  tan  3 **~ 


/ 


/ 


/ 


/ 


Zone  II 


/ 


Y 


Blade  2 


x 


Blade  1  '  * 

Figure  H-l.  Method  of  Singularities  Configuration 

The  solution  for  supersonic  flow  past  two  oscillating  blades 
in  Zone  II,  (Figure  H-l),  can  be  obtained  by  superposition  of 
single  blade  solutions.   For  this  purpose  the  disturbed   flow 
field  in  Zone  II  is  written  as  the  linear  superposition  of  the 
single  blade  solutions  for  blades  1  and  2.   However,  in  order 
to  satisfy  the  flow  tangency  condition  on  the  upper  surface 
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of  blade  2,  a  third  single  blade  solution  with  unknown  normal 
velocity  must  be  added.   Hence,  one  can  write 

*T  =  *x  +  *2  +  *Ip  (H_2) 

where 


-1  ,x-y/Mz-l  ,  .  -ik  M2  ,   r,     .„ 

'r^'«  v(c)e     »P=r(,-{)J0[,^(«-e)^^-l)y^de 

(H-3) 

y — - kM2 

/Mz-1  °  " 

(H-4) 

.  kM2 


*   =  _ 

IP  ^ 


—fQ11  va^)e      M2-l  1         J0[M^T/(xl-?)    (M  -1}yi  3dC 


(H-5) 


where  v   denotes  the  unknown  normal  velocity.   Equation  (H-3) 
represents  the  disturbance  due  to  the  first  blade,  equation 
(H-4)  the  disturbance  due  to  the  second  blade  alone,  and 
equation  (H-5)  the  interference  disturbance.   Note  that  the 
blade  2  and  interference  potentials  are  expressed  in  terms  of 
an  x, ,  y.  coordinate  system  with  origin  at  the  leading  edge 
of  blade  2 . 

Satisfying  the  boundary  condition  on  the  upper  surface  of 
blade  2  gives  an  equation  for  the  unknown  velocity,  v  ,  i.e., 

v  (x,y)=  —  /  *     v-U)e  M7-! 

°  /M^T  8y   ° 


•  J0[j^rr^(x-C)^-(M^-l)y^]d5 


at  y  =  d. 


(H-6) 
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Now  setting 

v(C)  =  ik  ZU)  +  Z'(£)  (H-7) 

and  assuming  the  blades  to  pitch  about  their  leading  edges, 
i.e.  , 

Z(C)  -  g  (H-8) 

and  then  performing  the  differentiation,  one  obtains 

/     \       1    I  ,X-y/Mz-l  r  .,„,>-\.».*/,-\n   ~i— 7 — 7-(x-£  ) 

v  (x,y)=      <(f      ■*      [lkZ  (€)+Z  (?)  ]e  Mz-1    s/ 


o 


/M^=T  1  o 


•       .  RMy      *  jt  ^m^T  '(x-5)*-(M*-l)y*]cU 
/(x-C)z-(Mz-l)yz     X  M   1 

_____         _____      _ .  kM2 

-  [ikZ (x-y/M/-l)+Z^ (x-y/Mz-l) ]  •  e"1^^^!7 

at  y  -  d  (H-9) 

Since,  in  general 

C_  =  ikcf>  +  4>  (H-10) 

P  x 

the  total  pressure  coefficient  can  be  written  as 

CP  ■  CP1  +  CP2  "  CPIp  <H-U) 

where  C  -,  is  Obtained  from  the  single  blade  solution  for 
blade  1  at  y  =  d, 

CD1=  JiL(/x-d/^^[-Z(U4z^(de-iM^r(x-° 

•  JqIj^t  /(x-OiMMM)da]dc4//J-d,^lr[-z'(5)4z''^)] 

kM2 
.e"aM2=T(x"5)  •  J0[M^I/(x'C)/'(MM)d'ldC 

+  ie-W  *  *  Vm^I  'x<-<M<-l)d<]]>  )         (H-12) 
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C  7    from  the  single  blade  solution  for  blade  2  at  y,  = 


k2     X  i  -j   kM  2  f..  -F) 


o, 


'p2   /m^T  o 


•  VM^T(xrU]d^  +  £  e"1*^  Xl   jo[m^i  xi][  }        (H~13) 

and 

cpIP  ■  ik  *ip  +  ♦«  (H-14) 

where 

i  •  kM2  /      a       i  u 

♦iP  "  j=£=   V  va(Y)e   iP=T        J0[FFTr(x1-Y)]dY 

(H-15) 


and 

a         -1    rxi    /  \  r  •  kM2   -j— -— — (x.  -£) 


.  kM 

*  J0[M^T(xl"Y)]dY     /01   a(Y)  e   FF"1   1     *  M^T 

•  J1[j^I(x1-Y)]dY?>  (H-16) 

where 

f  v—  kM2 

v  (x)  =   _l_;/>=-d'/M7::l[ikZ(U+z.(C)]e-iSrr[(x-C) 

/(x-C)2-(M2-l)dz     X  M  _1 

kM2 
-[ikZCx-d/M2^!)  +  z'  (x-d»/FF=T)]e":LM7^r  ( 

Note  that  the  substitution  of  v  (x)  into  equations  (H-15)  and 

(H-16)  for  <t>    and  4>     requires  one  to  relate  x  and  y   by 
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x  =  y  +  d  tan  B  (H-18) 

(see  Figure  H-l) 

These  integrals  were  computed  using  the  Trapexoidal  Rule 
and  IBM  subroutine  BESJ.   A  listing  of  the  program  follows 
in  Appendix  I. 
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APPENDIX  I 


PROGRAM  LISTINGS 


OOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOO 

•^  (nj  rn  ^  ir>  >o  p^  co  ^  o  ^  <nj  ro  vt  u^  so  r^  co  c>  o>-'r\irt'i>tu^vOr~- ooc^o--^  cxi  ro<r  mor-cocro— <<\jr<~)<MnvOr~oo 
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